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SUMMARY 

A new computational fluid dynamic code has been developed for the study 
of mixing and chemical reactions in the flow fields of ramjets and scramjets. 
The code employs an implicit finite volume, lower-upper symmetric successive 
overrelaxation scheme (LU-SSOR) for solving the complete two-dimensional 
Navier-Stokes equations and species transport equations in a fully-coupled and 
very efficient manner. The combustion processes are modeled by an 8-species, 
14-step finite-rate chemistry model whereas turbulence is simulated by a 
Baldwin-Lomax algebraic model. The validity of the code is demonstrated by 
comparing the numerical calculations with both experimental data and previous 
calculations of a cold flow helium injection into a straight channel and pre- 
mixed hydrogen-air reacting flows in a ramped duct. The code is then used to 
calculate the mixing and chemical reactions of a hydrogen jet transversely 
injected into a supersonic air stream. Results are presented that describe 
the flow field, the recirculation regions in front and behind the injector, 
and the chemical reactions. 


INTRODUCTION 

The recent interest in hypersonic vehicles has created a need for effi- 
cient numerical methods for the prediction of turbulent mixing and combustion 
in transonic and supersonic flows. In developing such a predictive method, 
numerical stability and efficiency are two important aspects that have to be 
considered. The stability problem arises mainly from the stiffness of the 
chemical source terms in the species concentration equations. The concern 
over the numerical efficiency, on the other hand, is mainly caused by the large 
number of species equations that have to be solved along with the flow equa- 
tions and the close coupling between the flow and the chemistry. In recent 
years a number of numerical methods (refs. 1 to 4), have been developed for the 
calculation of steady-state supersonic reacting flows. Most of these methods 
have treated chemical source terms implicitly to remove the stiffness associ- 
ated with the species equations. While this approach (treating chemical source 
terms implicitly) has been successful in circumventing the stiffness limita- 
tion, most existing methods have suffered from numerical inefficiency. The 
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majority of the reacting flow codes developed to date has adopted explicit 
integration methods in which the flow and the species equations are solved 
sequentially and are not directly coupled. Since in most chemically reacting 
flows the coupling between species and flow variables is strong, explicit meth- 
ods are generally not robust and very slow in convergence. Most implicit meth- 
ods, such as the popular Beam and Warming scheme (ref. 5), on the other hand, 
require the inversion of large block matrices. This can be very costly when a 
large number of species equations are solved in a coupled manner with the flow 
equations . 

In the present paper a new numerical method is described for the study of 
mixing and chemical reactions in transonic and supersonic flows. This method 
employs an implicit finite volume, lower-upper symmetric overrelaxation 
(LU-SSOR) scheme for solving the complete two-dimensional Navier-Stokes equa- 
tions and species transport equations in a ful ly-impl i ci t and fully-coupled 
manner, with real gas properties. Despite being implicit, the LU scheme 
requires only scalar diagonal inversions while most other implicit schemes 
require block matrix inversions. The use of scalar diagonal inversions offers 
large savings in computer CPU time. The LU scheme has recently been developed 
for nonreacting flows by the second author (refs. 6 and 7). Extensive tests 
have shown that this scheme is very robust and efficient for transonic and 
supersonic nonreacting flows. The objective of the present paper is to demon- 
strate the capability of the new LU scheme for chemically reacting flows. 

Several mixing and/or combustion problems are considered in this study, 
including a cold flow helium injection, premixed hydrogen and air reacting 
flows in a ramped duct, and a transverse injection of hydrogen into a super- 
sonic airstream with combustion. The chemical reactions are represented by an 
8-species, 14-step chemistry model of Westbrook (ref. 8) and Hitch et al . 

(ref. 9). The governing equations are Reynolds-averaged, and the closure of 
the equations is provided by the algebraic eddy viscosity model of Baldwin and 
Lomax (ref. 10). 


NOMENCLATURE 

Cp. constant pressure specific heat of species i 
C v . constant volume specific heat of species i 
D^j binary diffusion coefficient 
e total internal energy 

hi enthalpy of species i 

h? heat of formulation of species i 

' i 

I identity matrix 
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Mi molecular weight of species i 

N s number of species 

Pr^. turbulent Prandtl number 

p pressure 

q heat flux 

R u universal gas constant 

SCj. turbulent Schmidt number 

T temperature 
t time 

u streamwise velocity 

Ui streamwise diffusion velocity of species i 

v transverse velocity 

Vi transverse diffusion velocity of species i 

Xi molar fraction of species i 
x,y,z Cartesian coordinates 
Yi mass fraction of species i 
p laminar viscosity 
p denisty 

x stress tensor 
Subscripts : 

i ,j x and y grid indices 
Superscripts : 
n time index 

ANALYSIS 

Governing equations 

The two-dimensional unsteady Navier-Stokes and species transport equa- 
tions for chemically reacting flows, in Cartesian cordinates, are qiven bv 
(refs. 11 and 12): 
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( 1 ) 


where 


and 



0 

T xy 

T yy 

ux X y + VT yy -c * 
- pviY 

i 

'0 ‘ 

0 

S = 0 , 1 - 1,2, .... N s -1 

0 

Hi 

T xx ■ 2 t* u x - I f <u x ♦ v y> 

\y ~ f <u y * v x > 

T yy * 2 * iV y ' ! * <U x + U y> 
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N. 


P x = 


8T 


- k fi + p h i Y i u i 


q y = 


N s 

- k Vy + P i f 1 h i Y i *i 


The diffusion velocities are found by F i ck ’ s law 


Vi - -°im k Y 1 


Vi ' -D 1m k Y 1 


( 2 ) 


where D^ m = (1 - X. )/A (X./D^) is the effective binary diffusivity of 

species i in the gas mixture (ref. 12). The temperature and pressure are 
calculated iteratively from the following equations 


e = 


E Y h. 

1.1 ' ' 


^ + \ (u 2 + v 2 ) 
p 2 


(3) 


rT 


h. = h^ + 


C dT 
p i 


(4) 


'Ref 


N 


5 Y i 
P = pRuT ifi ^ 


(5) 


Thermodynamics and Transport Models 

The specific heat, thermal conductivity, and viscosity for each species 
are determined by fourth-order polynomials of temperature. The coefficients 
of these polynomials are supplied by McBride (ref. 13) and are valid up to 
temperature of 6000 K. The specific heat of the gas mixture is obtained by 
mass concentration weighting of each species. The thermal conductivity and 
viscosity of the mixture, however, are calculated using Wilke's mixing rule 
(ref. 14). 

The binary mass diffusivity j between species i and j is obtained 
using the Chapman-Enskog theory in conjunction with the Lennard-Jones intermo- 
lecular potential functions (ref. 14). The diffusion of a species in the gas 
mixture is approximated by Fick's law, i.e., treating the species i and the 
surrounding gas as a binary gas mixture, and the diffusion velocities for each 
species are calculated using equation 2. It should be noted that Fick's law 
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is a convenient approximation of the multicomponent diffusion equations 
obtained from the kinetic theory of gases, and the error introduced by the 
approximation may become significant if the pressure gradients in the flow 
field are very large or molecular weights of the species are very different 
(ref. 12). 


Chemistry and Turbulence Models 

In the present study, an 8-species, 14-step chemistry model is adopted 
for hydrogen-air reactions. This model is a reduced H 2 -air reaction scheme 
developed from a sensitivity analysis by Hitch et al . (ref. 9) performed on a 
more complete model originally proposed by Westbrook (ref. 8). This reduced 
reaction mechanism was shown to yield good agreement with the more complete 
model for both the ignition delay and the reaction times over a wide range of 
conditions (ref. 9). It was also shown in reference 9 that the 2-step global 
model developed by Rogers and Chinitz (ref. 15) predicted much shorter igni- 
tion times and longer reaction times in comparison to the 14-step model used 
in this report. 

The closure of the governing equations is provided by the Baldwin-Lomax 
(B-L) (ref. 10) algebraic eddy viscosity model and constant turbulent Prandtl 
and Schmidt numbers (Prt = Set = 0.9). In the formulation of the B-L eddy vis- 
cosity model, the distribution of vorticity is used to determine the length 
scales. The primary advantage of this model is that the boundary layer thick- 
ness which is often difficult to define in complex flows does not need to be 
calculated. However, it is found in the present study that in the vicinity of 
the injector very large vorticity is generated by the injector flow and, as a 
result, very large turbulent viscosities (about three orders of magnitude 
greater than the laminar viscosity) are calculated in this region. This large 
turbulent viscosity might not be physical and might affect the accuracy in pre- 
dictions of the flow structure as well as the mixing between the fuel and air- 
stream in the immediate vicinity of the injector. 

Due to the highly nonlinear nature of the chemical kinetics, effects of 
turbulence on reaction rates can be large, especially in regions with high tur- 
bulence intensity, such as the mixing layer, recirculation zone, and in the 
vicinity of the fuel injector. To account for such interaction effects would 
require a very sophisticated turbulence reaction closure or a direct numerical 
simulation (DNS). Since no effective turbulent combustion model is currently 
available which can handle detailed finite rate chemistry and the DNS methods 
for complex reacting flows are still under development, the interactions 
between turbulence and chemistry is not considered in the present study. 


NUMERICAL METHOD 
LU-SSOR Scheme 

Various numerical techniques have been used to solve the set of equations 
governing chemically reacting flows. Among these techniques, explicit schemes 
are generally slow in convergence when the flow involves high rates of chemi- 
cal reactions and heat release or large zones of recirculation. Most implicit 
schemes, on the other hand, require the inversion of banded block matrices and 
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become exceedingly expensive when the chemical system involves a large number 
of species. In the present study, the lower-upper symmetric successive overre- 
lation (LU-SSOR) scheme of Yoon and Jameson (refs. 6 and 7) is adopted to solve 
the coupled two-dimensional Navier-Stokes and species transport equations. The 
LU-SSOR scheme employs an implicit Newton iteration technique to solve the 
finite-volume approximation of the steady-state version of the governing equa- 
tions. Although the system of equations are formulated in a ful ly-impl icit, 
fully-coupled manner, the LU-SSOR scheme, unlike other implicit schemes, 
requires only scalar diagonal inversion for the flow equations and diagonal 
block inversion for the species equations. As a result, the LU-SSOR scheme 
has the advantage of a fast convergence rate while requiring a similar opera- 
tional count as that of an explicit scheme and hence is particularly attrac- 
tive for reacting flows with large chemical systems. The convergence of the 
Newton integration method is assured by the diagonal dominance of the coeffi- 
cient matrices of the LU-SSOR scheme. 

In the following, for simplicity, the derivation of the LU-SSOR scheme 
will be presented for the Euler equations. The final formulation for the 
Navier-Stokes equations will be given at the end of the derivation. A 
prototype implicit scheme for a system of nonlinear hyperbolic equations such 
as the Euler equations can be formulated as 


Q n+1 = Q n - (5 At |D x F ( Q n+1 ) + DyG (q" + 1 ) - S n+1 } 

-<1 - 13) At {d x F (o n j + D y G (q°) - S n } (6) 

where D x and Dy are difference operators that approximate 8/3x and 3/3y, 
and |3 is a positive number between 0 and 1. Here n denotes the time level. 
In this form, the scheme is too expensive since it calls for the solution of 

coupled nonlinear equations at each time step. Let the Jacobian matrices be 

A . | , B . | . and H . | 

and the correction be 

5Q = Q n+1 - Q n 

The scheme can be linearized by setting 

F(Q n+1 ) = F(Q n ) + ASQ + 0<||SQ|| 2 ) 

G(Q n+1 ) = G(Q n ) + BSQ + 0(||SQ|| 2 ) 

S(Q n+1 = S(Q n ) + HSQ + 0(||6Q|| 2 ) 

and dropping terms of the second and higher order. This yields 

{I + |3 At (D X A + DyB - H)} 5Q + At R = 0 (8) 
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where R i s the residual , 


R = D X F (Q n ) + D y G (Q n ) - S n 

If & = 1/2, the scheme remains second-order accurate in time, while for other 
values of [5, the time accuracy drops to first order. 

The unfactored implicit scheme (Eq. (8)) produces a large block banded 
matrix, which is very costly to invert and requires huge storage. Most resear- 
chers have adopted the ADI scheme which replaces the operator of equation (8) 
by the product of two one-dimensional operators 

Lx Ly SQ = -At R (9) 

This scheme requires relatively expensive tridiagonal or pentadiagonal block 
matrix inversions. Linear stability analysis has shown that ADI schemes are 
unconditionally stable for two-dimensional flows. In the three-dimensional 
case, however, the ADI scheme is unconditionally unstable. Although artifi- 
cial dissipation has some stabilizing effect, a large amount of dissipation 
can impair the accuracy. 

Jameson and Turkel (ref. 16) proposed the idea of a lower-upper (LU) fac- 
tored implicit scheme that is unconditionally stable in any number of space 
dimensions. There are many different ways to formulate an LU scheme (refs. 6 
and 7), and here we begin with the following form of equation (8) 


{I + |$ At (D A + + D*A + D B + + D*B - H>} SQ = -At R (10) 

x x y y 


where D and D are backward-difference operators and D* and D* are 

x y x y 

forward-difference operators. Here, two-point operators are used for steady 

flow calculations. A + , A”, B + , and B" are constructed so that the eigen- 
values of "+" matrices are nonnegative and those of matrices are nonposi- 
tive. The development of these matrices is extremely important for the suc- 
cess of LU-type schemes. One possibility, which yielded good results in the 
past (refs. 6 and 7) and is used in this work, is 


where 


A + = \ (A + v A I), A = j (A - v A I> 

B + = \ (B + v B I), B" = \ (B - v B I) 


( 11 ) 


v A > max ( | X A | ) , 


Vg > max ( |Xg| ) 


( 12 ) 


Here, X/\ and Xg represent eigenvalues of Jacobian matrices. After manipula- 
tion (the purpose of it will become apparent later), equation (10) can be 
written as 
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( 13 ) 


It 


I + 13 At (A + 


- A - ) + (B + 


B")] + 3 At [D x A + + D+ A + D y B + + D+ B" 
_ (A + - A") - (B + - B“) - h]U = -At R 


If we take "+" and matrices as given in equation (11), then, 


A+ - A- = v A I 
B + - B' = v B I 


Equation (13) can then be factorized to the following form 


j[l + 3 At (A + - A") + (B + - B')] + ^3 At <D X A + + D~B + - A + - B + - H)j| 

[l + 3 At (v A + v B )] _1 I |[l + 3 At (A + - A - ) '+ (B + - B')] 


+ [3 At (D+A - + D*B - + A + B ]j $Q = -At R (14) 


It should be mentioned that the factorization of equation (14) is not unique, 
and among various forms we have tried equation (14) gives the fastest conver- 
gence to steady-state solutions. If 3 = 1 , the scheme reduces to a Newton 
iteration in the limit At - 


(D~A + + D~B + - A' - B - - H) [(v ft + v^l]' 1 

(DV + D+B~ + A + + B + ) SQ = -R (15) 

x y 


It can be written as 

(D~A + + D~B + - A" - B" - H) (D*A~ + D+B~ + A + + B + ) SQ 
x y x y 

= -(v A + Vg) (D X F + D y G - S) (16) 


This equation can be further reduced to 
K * V B U - - H ] [ <V A 


For the Navier-Stokes equations, F and 
tion (17) are replaced by F - F v and 


+ V 1 + A Ul,j + B i,j + l] 60 

= -(v A + v B ) (D x F + D y G - S) (17) 

G on the right-hand side of equa- 
G - G v . That is. 
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(18) 


<V A + V 1 


- A t_l,j - B T,3-1 - H ] K * V B ) 1 * A Ul,j * 55 

- V - S ] 


= -(v A + Vg) [l 


D X (F - F v ) + D y (G 


Flux and Source Jacobians 

In the development of flux Jacobians A and B, evaluation of deriva- 
tives of pressure with respect to other independent variables is required. 

The exact formulation of these derivatives is quite complicated for flows 
involving real gases with nonequilibrium chemistry (ref. 17). An alternative 
which yields much simpler formulation is to use the concept of "equivalent 
Gamma." Let 



(19) 


Here y denotes the equivalent y, and h Sm and e Sm are the sensible 

enthalpy and internal energy (excluding the kinetic energy and heat of forma- 
tion) for the gas mixture. With the help of equation (19), the relationship 
for pressure can be written as 


p = (y - 1) pe, 


= (y - 1) 


pe - ph^ - £ 

m 


9 - £-(u 2 + / 


where h9 is the heat of formation of the gas mixture, 
r m 

N s 

h? = I Y. h9 

m i=l 1 r i 


(20) 


Since only N s - 1 species equations are actually solved, equation (20) needs to 
be rearranged to the following form 


P = (y - 1) 


pe - 


V 

l 

i = l 


P Y, h ° 



t - 2 


V 2 ) 


= <Y 


1) pe 



V 1 

E pY. Ah° - f-(u 2 + v 2 ) 
1-1 1 i 


( 21 ) 
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where 


Ah° « hj - h° 
f i f i f N s 

For real gases it is obvious that y = y(e,p,y^). However, as discussed by 
Golella and Glaz (ref. 18), this function of y varies very slowly with the 
change of other thermodynamic and flow properties, and hence we have treated 
y as a local constant in the evaluation of pressure derivatives. This approx- 
imation will not affect the accuracy of the converged steady-state solution. 

The final form of the Jacobian A is given as 


(Y - l)(f- - h°)-u 2 
* N 


(3 - Y)u 


-CY - 1 )v <Y - 1) -<Y - 1 )^h; -(Y - l) A h° 

r 1 r 2 


^ (Pe + p) + (v - 1 )ucf- - l (Pe + p) - (Y - l)u 2 -<Y - l)uv Yu 

p ^ N 


-uY } Y 1 00 

-uY 9 Y 9 00 


-uY m 


'Vi 


Vi 


-( Y - 1)u *h° -(Y - l)u 

' 1 T 2 


0 . 


-<Y - D A h? 


s-1 


-(Y - Du *t\° f 


( 22 ) 


where q2 = u^V. The form of B is similar to that of A. 
Equation (18) can be rewritten in the following form, 


[ DIAG i.j 




+ V 1 ♦ 



- -<v A * v B > [D X <F - F v ) + D <G - G v ) - S] (23) 

where DIAGi > j = [(v/\ + vg)I - H]jj are the (N s + 3) x (N s + 3) blocks in 
the diagonal of the matrix operator. These diagonal blocks are of the follow- 
ing form. 
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0 

0 


0 

0 


DIAG 


(V A * V 


9S 5 

*r 


9 s, 


N s+ 3 

~Zp— 


< V A + V 


aS 5 

3pu 


8 S 


V 3 


" 3pu~" 


(V A + V B ) 


3S 5 

ZpV 


< V A + V 


3S 5 

^Pe 




0 

0 

0 

0 

9 s c 


Bpv 


V 1 



<V A + V 


(24) 


In equation (24), 9S,/3p = 0 (since pyj are solved as independent varia- 
bles), and 3Sj /9(pu) and 3Si/3(pv) are taken as zero since they are usually 
small, although they can be evaluated via chain-rule, e.g., 


3S. /9SA fdj\ f 3e ^ 

3pU V 8T JpY . V3eJp,pY. l 8pu Jp,pv,pe,pY. 


where T denotes temperature, e = e 



n 



C v dT, 


and 


C = E Y.C . The terms of 3S,/3(pe) in equation (24) are evaluated using 
v m i=l v i 

a chain-rule formula similar to equation (25). The 3S i / 3( pYj ) terms can be 
easily calculated from the law of mass action (ref. 12). 


It is interesting to note that for nonreacting flows (S = 0 and H = 0 in 
eq. (23)) the present numerical method eliminates the need for banded block 
matrix inversions without using the diagonal ization procedure. In fact, with 
forward (for the first operator in eq. (23)) and backward (for the second oper- 
ator) sweeps in the diagonal directions (indices i and j increasing or 
decreasing simultaneously), only scalar diagonal inversions are needed to solve 
equation (23) for nonreacting flow problems. For reacting flows, due to the 
presence of the chemical source Jacobian H, the first operator on the left- 
hand side of equation (23) now requires block diagonal inversions. However, 
since in the present formulation the flow equations (continuity, momentum, and 
energy equations) have no source terms, the first four rows of the diagonal 
block (DIAG) of this operator have nonzero terms only in the diagonal. As a 
result, the first operator of equation (23) can be inverted in essentially two 
separate steps: the scalar diagonal inversion for the flow equations and the 

block diagonal inversion for the species equations. 
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Boundary and Initial Conditions 


For the problems to be considered here the inflows are always supersonic, 
so the upstream boundary conditions are provided by specifying the velocities, 
static pressure and temperature, and species mass fractions. For the super- 
sonic outflow the dependent variables are extrapolated from the interior. At 
the outflow plane, for grid points inside the wall boundary layers where the 
flow is subsonic, however, static pressure is specified (equal to the inlet 
pressure), and density is calculated using the ideal gas law. Along the solid 
wall, no-slip boundary conditions are specified. The wall is assumed to be 
adiabatic. The normal derivatives of pressure and species concentrations are 
also assumed to be zero. The injector is assumed choked (due to the high 
injection pressure), and hence at the injection slot all the dependent varia- 
bles are fixed at their initial values. Along a plane of symmetry, the normal 
derivatives of all the dependent variables are zero, except for the v veloc- 
ity, where v = 0 is specified. 

The governing equations are intialized by setting all the dependent varia- 
bles throughout the domain to the inflow conditions. 


RESULTS AND DISCUSSION 

One of the main purposes of this paper is to study the mixing and chemi- 
cal reaction of the transversely injected hydrogen in a supersonic airstream. 
This flow configuration is of particular interest to the design of scramjet 
engines because the transverse injection increases the fuel residence time in 
the combustor and the flow recirculation, resulted from the injection, can 
help to stabilize the flame. Since the present work is the first application 
of the LU-SSOR scheme to combustion problems, the ability of the method to com- 
pute complex reacting flows is assessed by first computing two simpler flows 
where either experimental data or previous numerical results are available for 
comparison. These two flows are: (1) the transverse injection of helium into 

a supersonic airstream, and (2) the combustion of premixed H 2 ~air supersonic 
flow in a ramped duct. 


Code Assessment 

An experiment was conducted by G.O. Kraemer and R.C. Rogers to study 
details of the flow field near a slot son leal 1 y injecting helium transversely 
into a ducted supersonic airstream. Results of this experiment were reported 
by Weidner and Drummond (ref. 19) along with their numerical predictions. The 
experimental apparatus is shown in figure 1. The conditions of the inflow air 
at the duct entrance are M = 2.9, T = 108 K, and p = 0.0663 MPa. The condi- 
tions of helium at the injection slot are M = 1.0, T = 217 K, and 
p * 1 .24 MPa. 

Results are given in figures 2 and 3. Figure 2(a) shows the static pres- 
sure contours. The helium jet partially blocks the cross flow resulting in a 
bow shock and strong pressure gradients around the Injector. The jet acceler- 
ates Immediately after injection because of the high injector pressure and 
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then forms several very strong jet shocks further downstream from the injec- 
tion slot. Figure 2(b) presents a comparison of the calculated and experimen- 
tal results for lower wall pressures upstream and downstream of the injection 
slot which is located at x = 0. The calculations slightly overestimate the 
peak pressure as well as the distance upstream from the slot injector where 
the pressure begins to rise. The pressure levels immediately downstream of 
the injector are also somewhat overpredicted. The cross-stream static pres- 
sures and helium mass fractions at 3.81 cm downstream of the injector slot are 
shown in figure 3. The peak pressure caused by the bow shock and the jet 
shocks is predicted to occur slightly closer to the lower wall than the experi- 
ment indicates. Overall, the agreement between the calculated results and the 
data are quite good. Considering the very complex flowfield associated with 
the transverse injection, the performance of the present numerical method is 
very encouraging. 

Next case considered for code validation is the combustion of the pre- 
mixed hydrogen-air supersonic flow in a ramped duct. The same flows have been 
computed by a number of CFD research groups (refs. 3, 4, and 20), and in the 
following results of the present method are compared with those of previous 
studies . 

The geometry and inflow conditions of the two test cases (T® = 900 K and 
1200 K) are illustrated in figure 4. The inflow is at temperatures below the 
ignition threshold. The viscous layer along the walls and the shock wave 
induced by the ramp increase the temperature to a value where significant reac- 
tions occur. The pressure contours for the case of T 0 = 1200 K are shown in 
figure 5(a). The leading edge shocks and the ramp shock are clearly seen. 

Also observed is the continual increase of pressure behind the ramp shock due 
to heat released during combustion. Comparisons of pressure and temperature 
between the present calculations and the results of references 4 and 20 are 
given in figures 5(b) and (c). The results are obtained along the y-station 
located approximately 0.13 cm from the lower wall. It can be seen that the 
agreement is reasonably good. The present calculation shows a higher level of 
temperature behind the ramp shock than indicated by the result of reference 20. 
This is attributed to the different chemical models used in the two studies. 

As discussed earlier, the reaction time of the 14-step H 2 -air chemistry model 
used in this study is shorter than the 2-step global model used in refer- 
ence 20, and, as a result, the present model predicts more complete chemical 
reactions and greater temperature rise behind the shock. Both the present 
results and the results of reference 20 show some oscillations in pressure 
immediately behind the ramp shock. This is caused by the differencing method, 
i.e., central differencing, used in both of the two studies. The TVD scheme 
used in reference 4 (fig. 5(c)) is seen to give a smoother pressure jump 
across the shock. The distributions of species mass fractions are compared 
next in figure 6. The present results show larger H 2 O and smaller OH concen- 
trations behind the shock than those of references 4 and 20, indicating faster 
and more complete reactions predicted by the present chemistry model. Again, 
this phenomenon is attributed to the difference in chemistry models used in 
these studies. The comparison illustrated in figures 2, 3, 5, and 6 indicates 
that the present results agree reasonably well with experimental data or previ- 
ous numerical calculations for flows considered here, and the LU-SSOR scheme 
has the ability to simulate complex flows found in scramjet engines. 
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Mixing and Combustion of H 2 Jet in a Supersonic Airstream 

To further demonstrate the capability of the present code, the mixing and 
combustion of a sonic transverse hydrogen jet injected from a slot into a Mach 
4 airstream in a two-dimensional duct is computed. The duct geometry, inflow 
conditions, and the main features of the flow structure are illustrated in fig- 
ure 7. The inflow conditions produce a global equivalence ratio of 0.71 
(hydrogen/air mass, ratio = 0.02). The convergence histories of the numerical 
calculations for this case and the other two cases discussed earlier are shown 
in figure 8. Also shown in this figure are the grid density and CPU time (on 
CRAY-XMP computer) for each case. Considering the fact that relatively compli- 
cated chemistry and thermodynamic and transport property models are used in 
these calculations, the convergence efficiency shown in figure 8 is certainly 
encouraging. 

The velocity vectors around the injector slot are illustrated in fig- 
ure 9. A large recirculation zone upstream of the injector can be clearly 
identified. The flow separation is caused by the adverse pressure gradient 
produced by the hydrogen jet and the shock-boundary layer interactions. There 
is also a smaller region of separation downstream of the injector caused by 
the blockage of flow by the hydrogen jet. In scramjet engines, the two recir- 
culation regions provide longer fuel residence times, better mixing of fuel, 
air, and hot combustion gas, resulting in better flame holding capabilities 
for the combustor. The static pressure contours are shown in figure 10. Both 
this figure and the schematic flowfield shown in figure 7 indicate that the 
transverse hydrogen jet partially blocks the axial flow and generates a strong 
bow shock just ahead of the injector. The hydrogen jet accelerates to super- 
sonic speed after injection and then decelerates by jet shocks not too far 
from the injector. The incidence and reflection (by the up-running bow shock 
from the bottom half of the duct) of the bow shock to and from the symmetric 
plane of the duct can also be seen in figure 10. The contour maps of the mass 
fractions of two major combustion products, i.e., H 2 O and OH, are presented in 
figure 11. The large mass fractions of H 2 O and OH found in the recirculation 
zone in front of the injector indicate that transverse injection of fuel pro- 
vides good flameholding capability for the simulated scramjet combustor. Fig- 
ure 11 also indicates that the injected fuel penetrates well into the main- 
stream, and good mixing and chemical reaction occur in the combustor. 


CONCLUDING REMARKS 

A computer code for analyzing two-dimensional chemically reacting flow 
fields has been developed. This code employs the LU-SS0R finite volume scheme 
which solves the Navier-Stokes equations and species transport equations in a 
fully implicit and coupled manner. The validity of this code is demonstrated 
by comparisons of the present calculations with experimental data and previous 
numerical results. The code has been used to simulate the flowfield in a 
scramjet combustor in which the hydrogen is transversely injected into a super- 
sonic airstream. The present results show that the new code can handle the 
complex flowfield associated with the fuel injection, mixing, and chemical 
reaction in a scramjet combustor. The new code has also been shown to be very 
efficient and robust for the flows considered in the present study. 
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Interactions between the turbulence and chemical reactions are not 
accounted for in the present study. The interaction effects are important in 
most ramjet/scramjet combustor flows, especially for simulating ignition, 
flameholding stability, and heat release rate, and should be considered in the 
future code improvement effort. 
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